Violation of the Einstein relation in Granular 
Fluids: the role of correlations 

A Puglisi^, A Baldassarri^ and A Vulpiani^ 

^ Dipartimcnto di Fisica, Universita La Sapienza, p.le Aldo Moro 2, 00185 Roma, 
Italy 

^ Dipartimento di Fisica and INFN, Universita La Sapienza, p.le Aldo Moro 2, 00185 
Roma, Italy 

E-mail: andrea.puglisiOronial . inf n. it 

Abstract. We study the linear response in different models of driven granular gases. 
In some situations, even if the the velocity statistics can be strongly non-Gaussian, 
we do not observe appreciable violations of the Einstein formula for diffusion versus 
mobility. The situation changes when strong correlations between velocities and 
density are present: in this case, although a form of fluctuation-dissipation relation 
holds, the differential velocity response of a particle and its velocity self-correlation 
are no more proportional. This happens at high densities and strong inelasticities, but 
still in the fluid-like (and ergodic) regime. 
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1. Introduction 

The transport properties of flowing dilute granular materials constitute an open problem 
in non-equilibrium statistical mechanics [HIS]. The existing kinetic theories aiming to 
deduce transport properties from the microscopic dynamics for the "usual" gases have 
a hard life, here, because of the presence of inelastic interactions among grains, that 
prevent the assumption of an equilibrium measure in the unperturbed state. Several 
approaches to this problem have been proposed in previous studies. 

In a large series of works, it has been considered a setup where energy is injected 
only through the boundaries of the system. In this treatment, the bulk is considered 
as a "freely evolving" inelastic gas, which would cool down if not driven by energy 
currents transported by the gas itself. It is possible to write down balance equations 
for local flelds such as density, velocity and kinetic temperature and, through a very 
delicate assumption of separation between microscopic and mesoscopic scales [3], the 
so-called granular hydrodynamics can be obtained |1] . Following these lines, fluctuation- 
response relations have been obtained with respect to the homogeneous cooling state, 
i.e. relaxation laws for small perturbations of a state whose fate is thermal death [5], [6]. 

A different approach consists in considering an alternative experimental setup where 
the starting state is much more similar to a thermal state, such as the thermodynamic 
equilibrium of a gas. Such a state can be prepared by coupling the energy source 
to all grains of the system, for example in granular materials fluidized by some air 
flow ^ , or otherwise in granular beds put on a vibrating plate |H1 E] . In such cases, the 
unperturbed fluid state is stationary and one can study how the system relaxes to it when 
a small perturbation is applied. Models for these granular stationary states have been 
proposed [lOl [HI [12] , showing the main differences with respect to the thermal state of 
a molecular gas: lack of equipartition, departure from Gaussian statistics of velocities, 
tendency to enhance spatial grain-grain correlations, clustering. It must be remarked 
that the obtained stationary state is intrinsically out of equilibrium: a net current of 
energy flows from the external source, through inelastic collisions, into heat [T3l [T4] . 
A Boltzmann equation for such class of granular fluids has been proposed [15], as 
well as a kinetic theory for transport coefficients [161 [13 • this theories have their 
validity in a suitable range of the parameters. Recently numerical studies have been 
performed showing that, in homogeneous situations, the Fluctuation- Response relation 
(FR) is valid in its near-equilibrium formulation, replacing the bath temperature with 
the internal granular temperature [18l[T9]. This has interesting consequences in the case 
of mixtures, where different components have different temperatures [20]: for instance, 
a linear response experiment on a massive tracer, performed to obtain a temperature 
measurement (a granular thermometer), yields the temperature of the tracer and not 
that of the surrounding gas. The verification of FR has been explained by means of a 
hydrodynamic approach by Garzo [21j, who connected it to the very small departures 
from the Maxwell-Boltzmann statistics. In the following, we discuss when this particular 
kind of Fluctuation Response relation ceases to be valid in a driven granular system: 
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it will appear that the most relevant ingredient is not the deviation from the Maxwell- 
Boltzmann statistics, but the degree of correlations among different degrees of freedom 
(d.o.f.), which increases as the total excluded volume decreases. 

The aim of this paper is to put this problem in the more general context of linear 
response theory for statistically stationary states, whose formulation has been given 
in |22l [23] and which can be described in very general terms. Consider a dynamical 
system X(0) = f/*X(0) whose time evolution can also be not completely 

deterministic {e.g. stochastic differential equations), with states X belonging to a 
A^-dimensional vector space. We assume a) the existence of an invariant probability 
distribution p(X), for which an "absolute continuity" condition is required, and b) the 
mixing character of the system (from which its ergodicity follows). These assumptions 
imply also that the system is time translation invariant (TTI). Now we introduce the 
two main ingredients of the theory: the response of the system to a small perturbation, 
and the time correlation of the unperturbed system that describes the relaxation of 
its spontaneous fluctuations. In the following we will indicate with (■) an average in 
the unperturbed system, i.e. weighting states with the invariant measure, and with 
(■) the time dependent average in the dynamical ensemble generated by the external 
perturbation. 

In the unperturbed system, the relaxation of spontaneous fluctuations is described 
by the time dependent cross correlations of two generic observables A(X) and -B(X) 

CMt) = (A(X(t))S(X(0))). (1) 

The average effect at time t on a variable Xi of a small external perturbation f{s), for 
instance on the variable Xj, applied at time s can be written, in the linear response 
regime, as: 



dX,{t) = / R,,,{t-s)f{s)ds, (2) 
Jo 

which defines the response function Ri j{t). The basic idea of the FR is to link the 
response functions {R} to suitable correlations {C} of the unperturbed system. 

For example, let us consider a colloidal particle in a fluid with friction constant 
7 and temperature T: one has a system with two d.o.f., position and velocity of the 
particle, {X, V), whose evolution is given by the Langevin equation: 

(3) 



dV' 2-fT 

where t] is a white noise, i.e. a Gaussian stochastic process with {ri{t)) = and 
{j]{t)ri{t')) = 6(t — t') (for the sake of simplicity we will assume the Boltzmann constant 
kB = 1). The time self-correlation of the particle velocity, Cvvit) = {V{t)V{0)) = 
(y^)e~'^*, when integrated from to cxd, determines the self-diffusion coefficient: D = 
{V{t)V{0))dt. It describes the asymptotic growth of the mean square displacement 
of the particle: {{X{t) — X(0))^/t) 2D. On the other side, when the momentum of 
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the particle is perturbed with a force f{t) = FQ(t), where 9(t) is the Heaviside step 
function, the response of the velocity itself, at very large times, reads: 

F r°° 

moo) = — / R{t)dt = fiF, (5) 



which defines the mobility yU. An easy computation gives 

/i = PD, (6) 

where (3 = 1/T. This relation, well known as the Einstein formula, obtained in his 
celebrated 1905 paper on Brownian motion [2U [25], is a primordial example of the 
Fluctuation Response relation: it relates, in fact, the response to a perturbation with 
the relaxation of spontaneous fluctuations. 

After the publication of the Einstein relation, a large amount of work [26 | [27| [28l [29| 
I5U] was devoted to generalize it to the class of (classical as well as quantal) Hamiltonian 
systems coupled to a thermostat at temperature T. This means considering dynamical 
systems with variables (q, p) whose time evolution is generated by a Hamiltonian Hq. 
The external perturbation appears as a perturbation of the Hamiltonian AT-C{t) = 
/(t)-B/(q, p). In this case, the Fluctuation Response relation can be written, among 
the others, in the following form: 
6A{t) 



^f^^^ RAjit) = f3{A{t)Bf{s)) . (7) 

The response of the observable A at time t with respect to a "force" / applied to the 
system at time s is related to the correlation, measured in the unperturbed system, 
between the observable itself at time t and the time derivative of the observable -B/ at 
time s. The latter is the one conjugated to / through the Hamiltonian. 

The fact that the FR theory was developed in the context of equilibrium statistical 
mechanics of Hamiltonian systems generated some confusion and misleading ideas on 
its validity. As a matter of fact it is possible to show that a generalized FR relation 
holds under the rather general hypotheses discussed above, i.e. basically the mixing 
property and the existence of an absolute continuous invariant measure p(X). The 
main result (for details see [23]) is the following Fluctuation- Response relation, valid 
when considering the perturbation at time of a coordinate Xj: 



From this relation all previous cases can be obtained. The Brownian motion of the 
colloidal particle, for example, has an invariant measure where position and velocity are 
independent: the part concerning V is of course a Gaussian with (y) = and (V^^) = j^. 
From formula (JHl), therefore, follows that 

Rvy = Mp{V{t)V{0)) (9) 

and this immediately returns the Einstein relation ([6]). In the rest of the paper, with a 
slight abuse of terminology, we will use the form "Einstein relation" to denote the time 
dependent form ([9]). 
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In the case of thermostatted Hamiltonian systems, on the other side, one has that 
p(q, p) oc exp(— /37i(q, p)). In such case equation ([8]) gives for example 

= /5(p^w|^^(o)) = -PJ^{p^iiM^)) ■ (10) 

In non Hamiltonian (and in general non-Gaussian) systems, the shape of p(x) is 
not known, therefore ([8]) does not give a straightforward information. However from it 
one can see that a FR relation still exists, stating the equivalence of the response to a 
suitable correlation function computed in the non perturbed systems. This mean that, 
from an ansatz on the invariant measure p one can directly deduce the response matrix. 

Following these lines, we analyze the response to small perturbations of a 
thermostatted granular gas, trying to connect the response properties of the stationary 
state with its many "anomalies" with respect to an equilibrium state. In particular we 
show (section 2) that, in homogeneous situations, even when the invariant measure 
is far from the Gaussian, the Einstein relation holds with good accuracy. On the 
contrary (section 3), when the granular effects (excluded volume and inelasticity) are 
strong enough to develop correlations between local density and velocities, the invariant 
measure of the system becomes highly non-trivial, and the Einstein relation is no more 
observed. We stress that the regimes considered here are always ergodic: this is a 
relevant difference with respect to previous studies on the violations of the Fluctuation- 
Response relation, which considered glassy systems in the non-ergodic (aging) phase [31] . 



2. Non Gaussian cases for which the Einstein relation holds 



2.1. The models 

We start to discuss the linear response of a dilute granular gas with grains of mass 
m = 1. Three different models, all in dimension d = 2, are considered here: 

(i) the homogeneously driven gas of inelastic hard disks in the dilute limit, evolving 
through stochastic Molecular Dynamics rules. 

(ii) the homogeneously driven gas of inelastic hard disks in the Molecular Chaos 
approximation, i.e. where its djTiamics is determined by Direct Simulation Monte 
Carlo (DSMC) algorithm 

(iii) the inelastic Maxwell model driven by a "Gaussian thermostat" 
In all the above models one has 

N d 

p{{^.,x,})=n^X[X[p^{vt^) (11) 

i=l a=l 

with n the spatial density n = N/V and Pv{v) the one-particle velocity component 
probability density function, v^"'^ the a-th component of the velocity of the i-th particle 
and d the system dimensionality. In particular, in models 2 and 3 this is true by 
assumption, while for model 1 it is well verified in simulations, as a consequence of being 
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diluted. In view of the fact that all discussed models are isotropic, in the following we 
will denote with v an arbitrary component of the velocity vector: the results do not 
change if v is the x or y component. 

The three models are known to display non-Gaussian Pv{v)- From the above 
discussion, it is expected that an instantaneous perturbation 6v{0) at time t = on 
a particle of the gas, will result in an average response of the form 



having defined Ci{t) = {v{t)v{0)) / (v^) 

Some previous studies already showed that for the inelastic hard disks model, in the 
dilute limit, it is very difficult to observe the discrepancy between R(t) and Ci(t), i.e. to 
see any "violation" of the Einstein relation for mobility and diffusion. In such studies, 
however, the deviation from a Gaussian Pv{v) was always small, i.e. consistent with a 
Sonine polynomial fit with a parameter 02 = j^pyi — 1^1 [15j. Actually, in a driven 
dilute system it is very rare to observe large departures from the Gaussian behavior. On 
the other side, studying models such as the thermostatted Maxwell model, or tuning the 
parameters of the DSMC algorithm for inelastic hard disks beyond the dilute limit, one 
can induce rather large deviations from the Gaussian, while condition (fTTl) still holds. 
Even if this may be far from being realistic, it is useful to assess the relevance of the 
non-Gaussian velocity pdf on the linear response of the gas. 

For the three models, the unperturbed dynamics is determined by a non-interacting 
streaming (where each particle is coupled to the thermostat only) plus inelastic collisions. 
For model (i), which is the closest to experiments of driven granular gases, the streaming 
part is described by the equations of motions of hard disks of diameter a = 1 moving 
in a square of area V = L ^ L with periodic boundary conditions and coupled to a 
thermal bath with viscosity 7 and temperature T;,: 



where a and a' indicate the Cartesian components. When two grains i and j touch, an 
instantaneous inelastic collision takes place, with a change of velocities given by 



where r G [0, 1] is the restitution coefficient (the elastic case corresponds to r = 1) and 
cr is the unit vector joining the centers of the two colliding particles. In the dilute limit, 
numerical simulations show that colliding particles are not correlated, i.e. Molecular 
Chaos holds. The system is known to display very different regimes when rf,/rc changes, 
where r?, = I/7 and = l/ciJc is the average mean free time between collisions of a 
single particle. If Tb/Tc ^ 1, the effect of collisions is very small and, even if inelastic, 




(12) 





(15) 



Violation of the Einstein relation in Granular Fluids: the role of correlations 



7 



the gas behaves as at equihbrium at temperature Tf,. In the opposite case, u/tc ^ 1, 
colhsions are dominant and the gas reaches a fluctuating stationary state with "granular 
temperature" Tg = -^^^^ < (^5 is smaller, the smaller r is). In this non-equilibrium 
regime, the velocity pdf is non-Gaussian with slow tails at very large \v\. Note that, 
here, Tc is not an external parameter, but is self-determined by the system: increasing 
na, i.e. reducing the mean free path, results in a smaller granular temperature and 
the same happens when increasing Th, therefore in both cases a direct increase of t^/tc 
is not obvious. A direct experience teaches that, when diluteness (volume fraction 
= nira'^/A <^ 1) is also required, then it is very difficult to obtain ^ r^,. For 
instance, when (p ~ 0.1, we usually observe Tc ~ O.lrf, or larger. 

Model (ii) consists of the same physical system but now the Molecular Chaos 
assumption is artificially enforced. This is achieved by disregarding the spatial 
coordinates of the particles, and selecting with a stochastic rule the pairs of particles 
involved in each collision. Time is discretized in small steps of length At (smaller 
than Tf, and Tc). At each step the discretized version of IHM is used to evolve the 
velocities of all the particles. Then, a number of collisions Nujc^t/2 is performed, 
where Uc = ^an^ixTg is the theoretical one-particle collision frequency for a dilute 
gas. Pairs i,j to collide are chosen with a probability proportional to the quantity 
— (vj — Vj) ■ crO(— (vj — Vj) • cr), with cr = (cos ^, sin 6') and Q chosen randomly with 
uniform probability in [0,27r). This process mimics the relative velocity dependent 
collision frequency in dilute gas of particles with hard core interactions. Since, typically, 
a particle i, after having collided with a particle j, will do a second collision with the 
same particle j after a number of collisions of order A^, any memory of the first collision 
will be lost and the new collision can be considered uncorrelated to the previous one. 
The phenomenology observed in this model is analogous to that of model (i) in its dilute 
limit. On the other side, here one can arbitrarily tune the ratio Tb/rc, increasing na and 
ignoring the inconsistency between the high density and the enforced Molecular Chaos. 

Finally, model (iii), the inelastic Maxwell model, is a further simplification of model 
(ii): couples of particles collide with a constant probabihty, i.e. independently of the 
relative velocity. The collision frequency is assumed to be = 1/A^. Moreover, the 
streaming part of the dynamics is performed with = and 7 = —A with A = (1— r^)/4, 
i.e. a negative friction and no random forces. This is the so-called Gaussian thermostat, 
which guarantees a constant kinetic energy, since in the non-driven case ( "free cooling" ) 
the total energy of the system would decay as ~ exp(— At). It has been shown that such 
"thermostat" is equivalent to consider the free cooling system and continuously rescaling 
all particles' velocity components by a factor The analysis of this idealized granular 

model is instructive for the following reasons. First, it has been shown [32] that such 
model has a stationary probability density function for the velocity hv with power law 
tails. In particular, it displays p„(f ) with high energy tails of the form v'^ with 6 = 4 in 
d = \ and 6 > 4 in ci = 2 (a good estimate for not too high inelasticity is h ~ 4/(1 — r)). 
Second, the simplification of the dynamics allows a direct analytical computation of 
time correlations and responses. 
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Figure 1. Left: correlation functions Ck{t) and response functions R{t) vs. time for 
different models described in the text. Right: pdfs of the x-component of the velocity 
(here denoted as v), for the same models. In all simulations N = 1000. In the MD 
simulation the box is of size 100 x 100, Tb = I/7 = 10 an = 1. All other parameters 
are reported in the figure. 



2.2. The numerical experiment 

The protocol used in our numerical experiment, for the three models, is the following: 

• First, the gas is prepared in a "thermal" state, with random velocity components 
extracted from a Gaussian with zero average and variance Tg{0). Positions of 
the particles, relevant only in model (i), are chosen uniformly random in the box, 
avoiding overlapping configurations. 

• Second, the system is let evolve until a statistically stationary state is reached, 
which is set as time 0: we verify that the total kinetic energy fluctuate around an 
average steady value and that this value does not depend upon initial conditions. 
In the case of model (iii) the energy is stationary by definition, therefore we ensure 
that a stationary velocity pdf is observed. 

• Third, a copy of the system is obtained, identical to the original but for one particle, 
whose X (for instance) velocity component is perturbed of an amount 6v{0). 

• Finally both systems are let evolve with the unperturbed dynamics. In models 
(i) and (ii), which involve random thermostats, the same noise is used. After a 
time tmax large enough to have lost memory of the configuration at time 0, a new 
copy is done with perturbing a new random particle and repeating the response 
measurement. 

This procedure is performed many times, in order to reduce the statistical errors 
for both the response and all required self-correlations in the unperturbed copy. In the 
following, averages indicated as (■) and (■) will have the meaning of averages over many 
realizations of this procedure. 

In Figured] we show the results of these experiments for the three different models. 
In the right frame the velocity pdfs are shown for different choices of the parameters in 
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those models. In Molecular Dynamics simulations of inelastic hard disks, even if quite 
inelastic, but still dilute (0 = mra'^/A < 0.1), the velocity pdf is not very far from a 
Gaussian, as in a DSMC with similar choices of the parameters (tc ~ O.lr;,). Increasing 
n in the DSMC leads to stationary regimes very far from thermal equilibrium, with 
Tg <C Tb and larger tails of the velocity pdf Pv{v), with 02 ~ 0.1. In view of relation (fT2ll . 
we have tried a three parameter fitting of the kind 



where cq is not independent because of normalization. In most of observed cases 
ka'f^maa;! |c2Wmaa:l5 with Vmax the largest value of V in the histogram. Therefore for our 
aims, in practice, we can drop the (negative) quartic term, retaining only ci and C2. The 
obtained fit appears to be very good, see Figure [H Using this formula in equation ( fT2l) . 
we get 



which also defines C2{t). This is an example of the discussion given in the introduction: 
an ansatz on the invariant measure leads to a link between the response function and 
some correlation functions. Note that here the ansatz is composed of two assumptions: 
positions and velocities of the grains are independent (this is exactly true for models 
(ii) and (iii)), equation (fTII) . and a specific non-Gaussian shape of Pv{v). Finally, for 
the Inelastic Maxwell Model the tail of the velocity pdf shows a power law decay with 
an exponent in agreement with its quasi-elastic limit 4/(1 — r) = 8 when r = 0.5. 

In the left frame of Figure [H, we have superimposed the response 6v(t)/6v{0) to the 
time self-correlations of different orders Ci{t) and C2{t) measured in the unperturbed 
system. From the above results, we learn that in all the considered models: 

• different correlations are almost identical (we do not show C^lt) = {v{t)v^{0)) / (f^), 
but the result is very close) 

• a very good agreement between R(t) and Ci(t) is observed, equivalent to a 
verification, within the limits of numerical precision, of the Einstein relation. 

The observation that the self-correlations Ck{t), at least for k = 1,2,3, are almost 
identical is very robust. With a precise statistics one can appreciate small differences 
at large times, proving that it is not an exact equivalence. Anyway, the measurement 
of response function is usually very noisy, and it is not easy to have a good signal/noise 
ratio at such late times. Therefore, for the practical purpose of the linear combination 
involved in the response, these small differences are negligible and the Einstein relation 
is practically satisfied. 

It is interesting to note that a rather similar situation is encountered when studying 
a different system, i.e. a gas of non-interacting particles whose velocities obey a Langevin 
equation with a non-quadratic potential: 



Pviv) = coexp(-cit;^ + C2\v^\ + csv^) 



(16) 



R{t) 



I = -2c,{vit)v{0)) + 3c2{v{t)\v{0)\v{0)) 
2c,{v^)C^{t) + 3c2{\vf)C2it), 



(17) 
(18) 




(19) 
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with U (f ) = Ci^;^— C2f |f P+C3f (with positive Ci, C2 and C3). A numerical inspection, not 
shown here, clearly indicates that Ci{t), C2{t) and C^lt) are almost indistinguishable. 

A simple condition can be given for the observed behavior. In fact, a generic time 
correlation for v{t) with a function /[f (0)] can be written as 



{v{t)f[v{0)]) = j dvtj dvoP,{vo)Vt{vt\vo)vtf{vo) = (20) 
dvop^{vo)f{vo){vt\vo), (21) 



where Vt{vt\vo) is the conditional probability of observing v(t) = Vt if f (0) = vq (time 
translation invariance is assumed) and {vt\vo) = J dvtVt{vt\vo)vt is the average of v{t) 
conditioned to v{0) = vq. 

If, for some reason, (ft|fo) = g{t)q{vo), with g and q two given functions, then the 
dependence on t results independent of the choice of the function /(f), i.e. on the order 
of the correlation. This happens in model (iii), where in spite of the non-Gaussian shape 
of the velocity pdf, the equivalence between R{t) and Ci(t) and any other correlation 

^ mm - m - e.p , (22, 



{i.(o)/wo)|) " "^V 4 

with any generic function / of the initial velocity value, is exact, see appendix B for the 
case d = 2 and [T^ for d = 1. 



3. Non homogeneous granular fluids 

The factorization of the invariant phase space measure, Eq. ( ITTl) . is no more obvious in 
model (i) when density increases. Correlations between different d.o.f., that is positions 
and velocities of the same or of different particles, appear also in homogeneously driven 
granular gases, as an effect of the inelastic collisions that act similarly to an attractive 
potential. Such a phenomenon has been discussed for this model of bath in [HI [33l [3l] 
and for other homogeneous thermostats in [TOl [121 [35]. [HI [23] it was also discussed 
the interplay between local density and local granular temperature, which in some very 
dissipative cases present strong fluctuations correlated to each other. These correlations 
indicate a breakdown of the factorization of the invariant measure, in particular at the 
level of velocity with respect to position of the same particle. As a matter of fact, 
these effects result in a strong violation of the Einstein relation, and in general of the 
equivalence between R{t) and Ci(t). 

Even in the presence of correlations, one can define and compute the marginal 
probability density function of the component x of the velocity of one particle i, 
projecting the phase space measure p({vj,Xi}): 

f^{v)= l[d^, n rfvfcrft;fp({v,x}). (23) 

k=l k=l,k^i 
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Figure 2. Left: correlations functions Ck{t) and response functions R{t) vs. time 
for a dense MD simulation: the response function is reported for different values of 
the perturbation Sv{0). Right: pdf of the x velocity component. The system has 
N = 1000, box of size 41 x 41, = I/7 = 10. In the simulation the mean free time 
between collisions is measured to be Tc = O.OSrft. 



However not necessarily this function has a role in the response function. For example, 
perturbing the x component of the velocity of the z-th particle and measuring the 
response of the same component, one obtains 

) ■ (24) 
t=o/ 

This is exactly what happens in model (i) when density is increased. In Figure[2], left 
frame, the correlation functions Ci{t) and C2{t) are shown, together with the response 
function measured with different values of the perturbation Sv{0). The very good 
agreement between different response functions guarantees that the system is indeed 
linearly perturbed. At the same time, the different correlations functions Ck{t) are very 
close, reproducing the phenomenology already observed in the previous dilute cases, 
with the difference that the time dependence is not exponential but slower, closer to a 
stretched exponential ~ exp(— (t/r)") with a < 1. Finally, looking at the velocity pdf 
of the gas, the previously proposed exponential of a cubic polynomial, Eq. ( |T6l) with 
a negligible C3 coefficient, is found to perfectly fit the numerical results. Therefore, if 
the correlations among the different d.o.f. are neglected, using equation ( fTTl) and the 
proportionality of the functions Ck{t), a verification of the Einstein formula R{t) ~ Ci{t) 
is still expected. The results displayed in Figure [21 left frame, demonstrate that this is 
not the case: the hypothesis of weak correlations among different d.o.f. must be dropped 
and the correct formula for the response is Eq. (^^. Unfortunately it is not very easy 
to use such a relation. 

The degree of violation of the Einstein formula increases with the volume fraction 
and the inelasticity 1 — r, as shown in Figure [3l where we have reported the ratio 
R{t) I C\ (t) as a function of time. This observation is consistent with the above argument: 



R{t) 



glnp({v,x}) 

dvf 



d\nf, 

dvf 



vt 
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Figure 3. Ratio between ttie response function R{t) and the normalized velocity self- 
correlation Ci{t). The ratio is 1 when the Einstein relation is satisfied. All the results 
come from MD simulations with TV = 1000 particles, Tb = 1 and Tb — 1/7= 10. 
Different values of the covered fraction and of the restitution coefficient are used, as 
shown in the figure. 



correlations among different d.o.f. increase when the probabihty of repeated contacts 
(the so-called "ring collisions" ) is enhanced, and this happens when the excluded volume 
and/or the post-collisional relative velocity are reduced. In the elastic case, r = 1, 
no violation is observed. A direct test of the existence of non-trivial correlations in 
the system is given in Figure |H Each point in this figure represents the value of 
CvN = (< vfNi > - < vf X Ni >)/(< vf X Ni >) for a given particle i, with 
Ni the number of particles in a squared box centered in Xj and of size L/lb, measured 
on a long trajectory of the unperturbed system. We observe that C^at increases together 
with the volume fraction 0. 



3.1. A Langevin model with two correlated variables 

In order to show in a clear way the role of correlations, we discuss now a simple model 
with only two variables: 

dx(t) 

miix(t) + mi2v{t) + aur]i{t) + cTi2T]2{t) (25) 



dt 

dv{t) 



mi2x{t) + m22v{t) + (T2irii{t) + cr22ri2{t) (26) 



dt 

If the matrices m and a are diagonal, the two variables are independent. Provided 
that the symmetric matrix m has negative eigenvalues and det (T 7^ 0, the pdf of (x, v) 
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Figure 4. Correlation between tlie square of tlie particle x velocity and local density 
CvN = (< vfN, > - < vf X Ni >)/(< vf >< N, >) for different values of the 
volume fraction, in MD simulations with N — 1000, = 1 and ti, = 10. In the 
dilute limit (jj —^ 0, C^n C'*^, which is different from zero because of the total finite 
number of particles. The red dashed line shows the estimate of C*^ obtained throwing 
N random velocities, extracted from a Gaussian distribution, into random boxes of 
the same size used in the MD, and repeating the measure over many independent 
realizations. 



relaxes toward a bi-variate Gaussian function. Instead of discussing the general form, 
we consider the case whose invariant joint pdf is 

2 2 

v) oc exp(-y - y + y)- (27) 

Of course the marginal pdf of each single variable is a Gaussian. Neglecting the 
correlation among x and f , the response of f to a perturbation on itself, would again be 
expected to be equal to Ci(t) = (f (t)f (0))/(f ^). On the contrary, the correct response 
is given using the full formula ([8]) applied to the joint pdf fl27|) . The result is 

R{t) = {v{t)vm-\mxm. (28) 

The difference between the Einstein formula and the correct response is shown in 
Figure [5] for a choice of the matrix rh. 



4. Conclusions 

In this paper we have reported the analysis of linear response in different models of driven 
granular gases, which have the property of rapidly reaching a statistically stationary 
state. The response function is directly related to the global invariant measure in 




Figure 5. Response R{t) and velocity correlation Ci{t) in the simple Langevin 
model with correlated variables discussed in Eq. (|25|) . with parameters mn = —1.1, 
mi2 = 0.8, TO22 = -1. 



the phase space, which is unknown for this kind of non-equihbrium systems. When 
positions and velocities of the particles are not correlated, the response of a perturbation 
on the velocity of a single particle is expected to depend on the singlet velocity pdf, 
which can be close or far from a Gaussian, depending on the model and on physical 
parameters. Nevertheless, the existence of a unique time scale that characterizes all 
possible correlation functions, makes the exact form of the velocity pdf irrelevant for 
the response function: the latter is, in practice, always indistinguishable from the 
normalized velocity self-correlation {v{t)v{0))/{v'^). This is equivalent to say that the 
Einstein relation is satisfied in all heated granular systems where correlations among 
particles are weak. On the other side, when excluded volume and energy dissipation 
occurring in collisions are increased, non-trivial correlations appear among positions and 
velocities of particles. The global invariant measure cannot be factorized anymore and 
the response function depends on it, i.e. on the specific parameters of the model. As a 
consequence, the Einstein relation is no more satisfied and the response function is not 
trivially predictable. In all simulations, the decay of R{t) is always faster than that of 
Ci{t): this is equivalent to state for the mobility that 

^ = _ / dtR{t) < -. (29) 
Jo 

Inspired by a recent work which reported violations of the Einstein relation in a 
non equilibrium model [36], we now conjecture an effective spatial dependence of the 
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Figure 6. Response R{t) and different correlation functions for the same MD 
simulation discussed in Figure O The normalized velocity self-correlations Ci {t) , as 
well as the correlation Cs{t) defined in Eq. ([3T|l . for different values of the coarse 
graining radius Lf,ox are reported. 



pdf of the velocity component for a particle at position x, at time t of the form 

p,(t;,x,t) ~exp|-L ^^r^j , (30) 

with m(x, t) a local velocity average, defined on a small cell of diameter Lhox centered 
in the particle. Such a hypothesis is motivated by the fact that, at high density or 
inelasticities, spatially structured velocity fluctuations appear in the system for some 
time, even in the presence of external noise fl2\ [35] . Following relation (JHl) we propose 
a formula for the response function to a velocity perturbation: 

R{t) =Cs = ^{vit){v{0) - u[xiO)]}). (31) 

Figure [6] shows that for small values of the coarse graining diameter Lbox (but still large 
enough to include 5 -i- 10 particles) relation fl?T]) is fairly verified. Note however that the 
proposed form (130!) cannot be exact, a spatial dependence of Tg should also be included. 
Furthermore, it is clear that ('u(x)) = for any point x, i.e. the local velocity field 
m(x) fluctuates in time. Thus, the above conjecture implies that the characteristic time 
of variation of these fluctuations is larger than the characteristic time of response of a 
particle: the particle feels, during its response dynamics, the "local equilibrium" average 
m(x). Further investigations are of course necessary to refine this promising argument. 
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5. Appendix A: the Inelastic Maxwell Model 

The dynamics of a particle in the Inelastic Maxwell Model is described by the following 
stochastic process: 

Av(t)At (with prob. 1 - At) 

Av(t) At + AvcoKv, u, cr) (with prob. At) ^ ^ 

where Avco/(v, u, cr) is the effect of an inelastic collision and has been defined in Eq. (JTSll . 

and u is the velocity of the collision partner. 

Let us define the two-times covariance matrix A^^{ti,t2) = {v^{ti)vy{t2)) ., with 
E {x,y}. Using the evolution law of the system, Eq. [321 one can calculate 



dt2 At2^o\ " At2 / 

= (^^^(^i) [^^-(^2) + ^v,oiA) = (34) 

= Kv^.{tl)v,{t2)) - ^^{v^{t^)a,[[w{t2) - u(t2)] ■ cr]) = (35) 

= Kvt.{ti)vu{t2)) (36) 

1 + r 

:^{Vf,{ti)[v^{t2)a^ - u^{t2)(y^ + Vy{t2)(yy - Uy{t2)ay\o^) = (37) 

1 + r 

= \{v^{ti)v^{t2)) —[{Vf,{ti)v^{t2)(T^a^) + {v^{ti)vy{t2)(Tya^) - (38) 

~{vi,(ti)U:r(TxCru) - {v^iti)uyaya„)] (39) 

which, assuming absence of correlation between pre-collisional velocities of different 
particles and between them and the impact vector cr, gives: 

1 + r 

MVf,{ti)v^{t2)) —[{Vf,{ti)v,x{t2)){axa^) + {v^{ti)vy{t2)){cXya^)] =(41) 

= f^A^u (42) 

where k = (A — ^-^) = — and we have used the fact that {a^a^) = \5^v Since 



the system is time translational invariant, we have that A^^{t) = A^[/(0) exp(K(t2 ~^i))- 
Now, one can perturb at a certain time the velocity v of a unique particle, with such 
a small perturbation that does not modify the rest of the system. Starting from (152]) . 
one easily computes the average response to such perturbation: 

= A(v> - i±^(F(v, u») = (a - i±^) (V). (43) 

where F(v, u, &) = {v^a^ - u^a^ + Vya^Oy - UyOya^, v^o^ay - u^a^Oy + Vya^ - Uya^) and 
we used the fact that (u) = 0. The result is that the average response decays as the 
velocity- velocity correlation. 



v(t+At)-v(t) 
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A more general result can be also obtained, starting from the stochastic evolution 
equation fl32l) . which can be rephrased as 

v(t + At) = i^(t)v(t) + J(t) (44) 

where Kit) and J(t) are two uncorrelated stochastic (two-dimensional) matrices, with 
{K^yit)) 7^ and (Jijit)) = . Starting at time with v(0) = Vq and iterating Eq. (jH]), 
one finds that (v(t)|vo) = I/*vo where L^^ = {K^jf). From relation (|20l) (generalized to 
two dimensions), the proportionality of all the correlation functions follows. 

Acknowledgments. - We wish to thank F. Cecconi, M. Cencini, L. Leuzzi and U. Marini 
Bettolo Marconi for discussions and a critical reading of the manuscript. 
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